#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <time.h>
#include <pthread.h>
#include <gsl/gsl_statistics_double.h>
#include "random.h"
#include "test-bench.h"

#define DEBUG 1
#define SAMPLE_SIZE 100000
#define THREAD_COUNT 128 /**< Number of threads to use in the multithreaded test. */

void random_setup_test() {
  setup_random_environment();
  setup_rng("exponential", cJSON_Parse("{\"type\":\"exponential\", \"lambda\":50}"));
  setup_rng("exponential1", cJSON_Parse("{\"type\":\"exponential\", \"lambda\":65}"));
  setup_rng("exponential2", cJSON_Parse("{\"type\":\"exponential\", \"lambda\":75}"));
  setup_rng("exponential3", cJSON_Parse("{\"type\":\"exponential\", \"lambda\":85}"));
  setup_rng("exponential4", cJSON_Parse("{\"type\":\"exponential\", \"lambda\":95}"));
  setup_rng("exponential5", cJSON_Parse("{\"type\":\"exponential\", \"lambda\":105}"));
  setup_rng("exponential6", cJSON_Parse("{\"type\":\"exponential\", \"lambda\":115}"));
  setup_rng("exponential7", cJSON_Parse("{\"type\":\"exponential\", \"lambda\":125}"));
  setup_rng("exponential8", cJSON_Parse("{\"type\":\"exponential\", \"lambda\":0.02}"));
  setup_rng("constant", cJSON_Parse("{\"type\":\"constant\", \"value\":100}"));
  setup_rng("uniform", cJSON_Parse("{\"type\":\"uniform\", \"min\":0, \"max\":100}"));
  setup_rng("normal", cJSON_Parse("{\"type\":\"normal\", \"mean\":100, \"std\":25}"));
}

int test_rng() {
  return 0;
}

int verify_rng(char *string, double target_mean, double target_std, double target_skew, double target_kurtosis) {
  double samples[SAMPLE_SIZE];
  double mean, std, skew, kurtosis = 0;
  int result = 0;
  int i = 0;
  for(i; i < SAMPLE_SIZE; i++) {
    samples[i] = get_random_number(string);
  }
  mean = gsl_stats_mean(samples,1,SAMPLE_SIZE);
  std = gsl_stats_sd_m(samples,1,SAMPLE_SIZE,mean);
  skew = gsl_stats_skew_m_sd(samples,1, SAMPLE_SIZE, mean, std);
  kurtosis = gsl_stats_kurtosis_m_sd(samples, 1, SAMPLE_SIZE, mean, std);
  //Check against the target values
  result =  test_matches_val(mean,target_mean,(target_mean*0.02)) && test_matches_val(std,target_std,(target_std*0.02));
  if(target_std != 0) { /* Check skew and kurtosis if we want a non zero standard deviation. */
    result = result && test_matches_val(skew, target_skew, 2) && test_matches_val(kurtosis, target_kurtosis, 2);
  }
  if(result != 1 && DEBUG) {
    printf("FAILURE:\n");
    printf("Mean:\nTarget:%6.2f%5s%6.2f\n", target_mean,"Actual:", mean);
    printf("Standard Deviation:\nTarget:%6.2f%5s%6.2f\n", target_std,"Actual:", std);
    printf("Skew:\nTarget:%6.2f%5s%6.2f\n", target_skew,"Actual:", skew);
    printf("Kurtosis:\nTarget:%6.2f%5s%6.2f\n", target_kurtosis,"Actual:",kurtosis);
  }
  return result;
}

int random_single_threaded_test() {
  return !(verify_rng("exponential",50,50,2,6)
           && verify_rng("exponential1",65,65,2,6)
           && verify_rng("exponential2",75,75,2,6)
           && verify_rng("exponential3",85,85,2,6)
           && verify_rng("exponential4",95,95,2,6)
           && verify_rng("exponential5",105,105,2,6)
           && verify_rng("exponential6",115,115,2,6)
           && verify_rng("exponential7",125,125,2,6)
           && verify_rng("exponential8",0.02,0.02,2,6)
           && verify_rng("constant",100,0,0,0)
           && verify_rng("normal",100,25,0,0)
           && verify_rng("uniform",50,28.86,0,-1.2));
}

void * random_thr_test(void *args) {
  return (void *) random_single_threaded_test();
}

int random_multi_threaded_test() {
  pthread_t threads[THREAD_COUNT];
  int i = 0;
  for(i; i < THREAD_COUNT; i++) {
    pthread_create(&threads[i], NULL, random_thr_test, NULL);
  }
  for(i=0; i < THREAD_COUNT; i++) {
    void *p;
    pthread_join(threads[i], &p);
    if((int)p) {
      return 1;
    }
  }
  return 0;
}

int main() {
  printf("Starting Random Number Generator Tests:\n");
  wg_test(random_setup_test,"Random Numbers Singled Threaded", random_single_threaded_test);
  wg_test(random_setup_test,"Random Numbers Multi Threaded", random_multi_threaded_test);
  return 0;
}
